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Excited states of the of the double-stranded DNA model (A)i2-(T)i2 were calculated in the frame- 
work of the exciton theory. The off-diagonal elements of the exciton matrix were calculated using 
the transition densities and ideal dipole approximation associated with the lowest energy tttt* excita- 
tions of the individual nucleobases obtained from TDDFT calculations. The values of the coupling 
calculated with the transition density cubes (TDC) and ideal-dipole approximation (IDA) meth- 
ods were found significantly different for the small inter-chromophore distances. It was shown that 
the IDA overestimates the coupling significantly. The effects of the structural fluctuations were 
incorporated by averaging the properties of the excited states over a large number of conformations 
obtained from the MD simulations. 



I. INTRODUCTION 

DNA a remarkable carrier of code of life, is very sta- 
ble with respect to the photochemical decay. The path 
chosen by Nature to protect DNA is through the very 
rapid decay pathways of the electronic excitation energy. 
Given the importance of DNA in biological systems and 
its emerging role as a scaffold and conduit for electronic 
transport in molecular electronic devices, |T] DNA in 
its many forms is a well studied and well characterized 
system. What remains poorly understood, however, is 
the role that base-pairing and base-stacking plays in the 
transport and migration of the initial excitation along 
the double helix. [2j|3l|4] 

The absorption of UV radiation by DNA initiate a 
number of photochemical reactions that can ultimately 
lead to carcinogenic mutations. [3 |6l [3 [3 [9] The 
UV absorption spectrum of DNA largely represents the 
weighted sum of the absorption spectra of it constituent 
bases. However, the distribution of the primary photo- 
chemical products of UV radiation, including bipyrim- 
idine dimers, [TD] is depends quite strongly upon base 
sequence, which implies some degree of coupling between 
the DNA bases. [3J Inasmuch as both the base stack- 
ing and base pairing are suspected to mediate the ex- 
cess of electronic excitation energy, understanding of the 
excited-state dynamics is of primary importance for de- 
termining how the local environment affects the forma- 
tion of DNA photolesions. 

Recent work by various groups has underscored the 
different roles that base-stacking and base-pairing play 
in mediating the fate of an electronic excitation in DNA. 
[5J[3] Over 40 years ago, Lowdin discussed proton tunnel- 
ing between bases as a excited state deactivation mech- 
anism in DNAjTT] and evidence of this was recently re- 
ported by Schultz et al. [12] In contrast, ultrafast flu- 
orescence of double helix poly (d A)- poly (dT) oligomers 
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by Crespo-Hernandez et oZ. [2] and by Markovitsi et al. 
[5] give compelling evidence that base-stacking rather 
than base-pairing largely determines the fate of an ex- 
cited state in DNA chains composed of adenosine and 
thymine bases with long-lived intrastrand states form- 
ing when ever adenosine is stacked with itself or with 
thymine. However, there is considerable debate regard- 
ing whether or not the dynamics can be explained via 
purely Frenkel exciton models [U |13j [14] or whether 
charge-transfer states play an intermediate role. [TS] 

Upon UV excitation, the majority of excited molecules 
shows a subpicosecund singlet lifetimes. [IH] [TTj [HI [H] 
Owing to the technical difficulties in measuring the ul- 
trashort lifetimes the study of the charge and excitation 
energy transfer in DNA has only recently received much 
of attention with the advances in the femtosecond spec- 
troscopy. Although, so far, no clear picture of the excited 
-state deactivation mechanism has been offered by the 
experiment, two possible decay channels have been in- 
vestigated. Kohler and coworkers in their recent study of 
the duplex poly (dA) -poly (dT) suggested that 7r-stacking 
of the DNA base determines the fate of a singlet elec- 
tronic excited state. Alternative decay mechanism in- 
volves interstrand hydrogen or proton transfer. Douhal 
and coworkers observed excited-state proton transfer in 
base pair mimincs in gas-phase. |20j The experimental re- 
sults suggests that these very fast decay pathways play an 
important role in quenching the reactive decay channels 
and providing DNA with intrinsic photochemical stabil- 
ity. However, they do not provide a clear picture which 
arrangement of bases, pairing or stacking, is of primary 
importance. 

Until recently, most theoretical investigations of exci- 
tation energy transfer in DNA helices has been within the 
Frenkel exciton model which treats the excitation as a co- 
herent hopping process between adjacent bases. [211 [22] 
This model has tremendous appeal since it allows one 
to construct the global excited states (i.e of the com- 
plete chain) in terms of linear combinations of local ex- 
cited states. The key parameter in the evaluation of the 
electronic excitation energy transfer (EET) is the elec- 
tronic coupling between the individual bases. To a first- 



2 



order approximation, the base to base coupling can be 
estimated using a dipole-dipole approximation in which 
the interaction between the donor and acceptor is cal- 
culated using only the transition dipole associated with 
each chromophore. While this approach is certainly suit- 
able for cases in which the distance between the donor 
and acceptor sites is substantially greater than the molec- 
ular length scale. In case of double stranded DNA, where 
the DNA bases are in relatively close contact compared 
to their dimensions this approach leads to the neglect of 
the effect of the size and spatial extent of the interacting 
transition densities associated with each chromophore. 

By far the most precise way to calculate the coupling 
elements is to directly integrate the Coulomb coupling 
matrix element between transition densities localized on 
the respective basis. [23^ The accuracy is then limited only 
by the numerical quadrature in integrating the matrix 
element and by the level and accuracy of the quantum 
chemical approach used to construct the transition den- 
sities in the first place. Futhermore, one must perform 
a quantum chemical evaluation of the coupling elements 
between each base at each snapshot along a molecular dy- 
namics simulation in order to properly take into account 
the fluctuations and gyrations of the chain itself. This 
is a formidable task, one that has prevented an accurate 
benchmarking of the excited state electronic structure of 
realistic DNA chains. 

In this paper, we present the results of simulations and 
calculations of accurate interbase exciton couplings for A- 
T strands of DNA in water in an attempt to provide such 
a benchmark. Starting from a molecular dynamics simu- 
lation of a model DNA sequence in water at the correct 
salt concentrations, we mapped out the evolution of the 
phot ochemic ally relevant excited states within a Frenkel 
exciton model in which the couplings were computed us- 
ing both the ideal dipole-dipole approximation (IDA) and 
using the transition density cube approach (TDC).[23] 



A. Exciton model 

The excited states of the {A) 12 ■ (T)i2 were calculated 
in the framework of the exciton theory j24l |25]. In this 
approach the total Hamiltonian for the super system of 
N molecules is written as the sum of N Hamiltonians of 
isolated molecules Hn and the intermolecular interaction 
potential Vnm between the molecules n and m. 
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The singly excited states of the system are described in 
term of N locally excited configurations 



(2) 



where 0; corresponds to the excited state wavefunction of 
the chromophore I whereas all the other molecules m are 
in their ground state 0„. denotes the corresponding 
wave function of the super system. Consquently, the ex- 
citon states of the supramolecular system can be written 
as a linear combinations of the excited states localized on 
each chromophore. 
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The diagonal elements of the exciton matrix 
{^l^\H\^l^) are simply excitation energies of chro- 
mophore n from its ground to i*'' excited state, Sq — > 
Si. The off-diagonal elements (^'Ji? j*!':^) written as 
('^n^ml^l'^n^m) corrcspoud to cxcitou couphng. It can 
be interpreted as the electrostatic interaction energy be- 
tween the transition densities corresponding to 5*0 —> Si 
and 5*0 — > 5*^. 

A measure of delocalization of the exciton states can 
be obtained from the inverse participation ratio (IPR) 
(1/Lfc) which represents the number of coherently cou- 
pled chromophores in a given eigenstate k. In the general 
case with more than one electronic transition per chro- 
mophore, Lfc is written as follows: 



II. METHODOLOGY 



The calculation procedure consisted of several steps. 
In the first stage the molecular dynamics (MD) calcula- 
tions were carried out to sample a range of conformations 
of (A) 12 -(T) 12 model of DNA double-helix. The transi- 
tion densities of the individual nucleobases obtained from 
time dependent density functional theory (TDDFT) cal- 
culations were subsequently superimposed with the in- 
stantaneous conformations from the MD simulations in 
order to calculate the coupling between the electronic 
transitions of the individual bases. In the final step, the 
excited-states of the model were calculated within the 
Frankel exciton model. 
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where k denotes a given eigenstate and i an electronic 
excited state of a chromophore. 

For purposes of developing a model, we can cast the ex- 
citon Hamiltonian as a SU{2) (g) SU{2) lattice model 
consisting of localized hopping interactions for exctions 
between adjacent base pairs along each strand (taj) as 
well as cross-strand terms linking paired bases {hi) and 
"diagonal" terms which account for the tt stacking inter- 
action between base j on one chain and base j ± 1 on 
the other chain {rf) in which rj denotes coupling in the 
5'-5' direction and coupling in the 3'-3' direction. Fig. 
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FIG. 1: Schematic view of the exciton couphng terms for 
Eq-i 



[T] shows a schematic view of the various coupHng terms 
between each nucleotide base. 

j 

where -i/jj and ipj are SU{2) spinors that act on the 
ground-state to create and annihilate excitations on the 
jth adenosine or thymidine base along the chain. The 7 
operators are the 2x2 Pauli spin matrices with ijjj = 71 ipj 
and 7+ -|- 7_ = 71 providing the mixing between the two 
chains. 

Taking the chain to homogeneous and infinite in ex- 
tent, one can easily determine the energy spectrum of 
the valence and conduction bands by diagonalizing 
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where sa.b ^md t^^r Sire local excitation energies and 
intra-strand hopping integrals, h is the coupling between 
Watson-Crick bases. When the interchain diagonal cou- 
plings are equal, r"*" = r~, Eq.|6]is identical to the Hamil- 
tonian used by Creutz and Horvath [27] to describe chi- 
ral symmetry in quantum chromodynamics in which the 
terms proportional to r are introduced to make the "dou- 
blers" at g cx TT heavier than the states at g oc since the 
off-diagonal coupling is now momentum dependent. 

One of the serious deficiencies with this model as it 
stands thus far is that for DNA each of the interactions 
described is very sensitive to the geometric fiuctuations 
of the DNA chain itself. [T3J [H] Hence, we need to con- 
sider each of the couplings as being parametrically depen- 
dent upon the instantaneous molecular geometry of both 
the individual bases and the chain itself. This is assum- 
ing there is no additional contribution from the solvent 
and ions surrounding the DNA chain. Assuming that 



the electronic time scale is fast compared to the typical 
time scale for geometric fluctuations of the DNA chain 
(lO^i-* - IQ-^^s for longitudinal and lO^^^ _ ;lO"12 g foj. 
the lateral motions of bases in DNA double helices[28 ), 
we can consider at least the initial electronic dynamics 
as occurring in a fixed nuclear framework and subsequent 
dynamics as adiabatically following the nuclear motion. 
Nonadiabatic contributions can not be completely dis- 
counted; however, the dominant non-adiabatic couplings 
are intermolecular in origin or involve proton between 
adjacent bases. [Tni^lHO] 



B. Transition densities and interactions 

1. Exciton- exciton interactions 

Each off-diagonal term in our Hamiltonian of Eq.[5]can 
be calculated according to 
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where the two terms in the numerator, and 
are the three dimensional charge distributions (transi- 
tion densities) associated with the ground and electronic 
excited states i and j of molecules n and m, respectively, 
with the separation between the elements a and b equal 
to rab- The corresponds to the electrostatic repul- 

sion energy between the two charge distributions AI^'' and 
of isolated chromophores. The calculations of the 
Coulombic couplings using the three dimensional charge 
distribution takes into account the size and the spatial ex- 
tent of the transition density and is valid at all molecular 
separations as opposed to the ideal dipole approximation 
(IDA). In the latter only the dipole moment of the transi- 
tion density is considered for calculations of the couphng 
terms which makes the computations of the off-diagonal 
elements much more efficient. However, this approxima- 
tion breaks down at the small donor-acceptor separations 
for which the spatial extent of the transition density be- 
comes important. 

To account for the dynamics of the DNA chain itself, 
we performed a series of molecular simulations of the 12 
base pair duplex DNA (AT) (Figurej2| with about 12,000 
water molecules and counter ions. |50j Once the system 
was minimized and equilibrated in the NVE ensemble at 
300K, we integrated the dynamics for an additional 80 
ps, sampling the DNA configuration every 10 fs. Even 
though we are dealing with a relatively small strand, it 
remains too large for an accurate evaluation of its elec- 
tronic structure. Consequently, we make the approxi- 
mation that the excited states of the molecule itself can 
be written as a linear combination of excited states local- 
ized on the instantaneous positions of each base along the 
chain. Furthermore, given the computational cost asso- 
ciated with evaluating the excited states of even a small 
molecule, it is prohibitive to perform such calculations 
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for each base at each time-step. Our approach, then, 
is to perform an accurate evaluation of the local transi- 
tion densities based upon the geometries of the isolated 
DNA bases, then map these densities onto the instanta- 
neous positions of the bases from the molecular dynamics 
simulations (Figure [2|. From this, we can evaluate the 
exciton-exciton coupling (Eq. [?]) in which Af"* and AfJ^^ 
are the transition densities about the instantaneous po- 
sitions of bases n and to. 



2. Excited states of individual bases 

The geometries of the DNA bases, adenine, guanine, 
cytosine, and thymine in their most common tautomeric 
forms were optimized at the MP2/TZVP level of the- 
ory in chloroform using GaussianOS suite of programs. |31j 
The optimized geometries were subsequently used to cal- 
culate the singlet excitation energies in gas phase at the 
TD-DFT level using PBEO functional and TZVP basis 
set augmented with the diffusion functions on all atoms 
as implemented in ORCA.[35] Additionally, the excita- 
tion energies were also calculated for the standard nucle- 
obase geometries obtained from the X3DNA.[331 In these 
calculations the deoxyribose and phosphate groups were 
replaced with hydrogens using the Chimera program. |34j 
Without further optimization of the structures, the exci- 
tation energies were calculated at the same level of the- 
ory as used before for the MP2 optimized structures. 
Fig. [3] shows both the transition density and direction 
of the transition dipole moment for each base as given 
by TDDFT after optimization at the MP2/TZVP level 
in a CHCI3. Transition moments were calculated using 
TDDFT with PBEO functional and aug-TZVP basis set 
in vacuum. The calculated excitation energies are sum- 
marized in Table HI 

The transition densities associated with the allowed 




FIG. 2: Schematic structure of the (A)i2-(T)i2 oligomer used 
in the MD simulations (left). The residue numbering is shown 
in the middle. The graphic on the right shows the three di- 
mensional representations of densities corresponding to lowest 
energy tttt* transitions of adenine and thymine superimposed 
with residues 7 and 19, respectively, of the (A)i2-(T)i2 model. 




FIG. 3: The transition densities and transition moments in 
the nucleobases: adenine, thymine, guanine, and cytosine. 
For purines the solid and dotted lines indicate the direction 
of the transition dipoles associated with the first and second 
lowest energy tttt* transitions, respectively. 

TTTT* excitations of the individual nucleobases, defined by 

M,g{r;s)^\^g){^Mr (8) 

where e and g correspond to the excited and ground 
states of the chromophore, were calculated using ORCA 
program. The densities are written in form of a charge 
distribution over three-dimensional grid of points, such 
that the integrated charge vanishes, according to 

rz+S^ py+Sy rx+S^ 

Meg (x, y,z):^Vs / / «'g^':ds dx dy dz{9) 

<y z J y J X 

where Vs = S^SyS^ is the element volume and the S^., Sy, 
Sz are the steps along the coordinate axes. The grid size 
has to be a compromise between the accuracy and the 
speed. Denser grids render the calculations very taxing 
while too small grids introduce large errors in the calcu- 
lated coupling elements. A satisfactory compromise was 
obtained for the cube files with 40 voxels along each axis 
{x, y, and z) which corresponds to total number of 64000 
elements. In case of the single nucleobase the volume of a 
single element is than 0.03 A'^. The changes in the mag- 
nitude of the coupling calculated for cubes with number 
of elements larger than 64000 was below 0.1 cm~^ . Ow- 
ing to the finite size of the cube the integrated charge 
over space was not exactly zero. The residual charge for 
all the transition density cubes was below O.Ole and was 
compensated by adding equal amount of charge to each 
volume element to bring the integrated charge over the 
cube volume to zero. 
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TABLE I: Vertical singlet excitation energies (eV) followed by oscilator strength in parentheses of the lowest electronic transi- 
tions of the nucleobases calculated using TDDFT at PBEO/aug-TZVP level of theory and MRCI in vacuum for the standard 
geometries taken from SDN A and optimized at MP2/TZVP level in chloroform. For the calculated excited state energies the 
transitions with a tttt* character are indicated with a boldface. 





Method 


Geometry 


So Si 


So —> S2 


So —^ S3 


So — > Si 


So — > S5 


A 


TDDFT 


MP2 


5.00 (0.002) 


5.29 (0.230) 


5.38 (0.069) 


5.49 (0.009) 








standard 


5.17 (0.001) 


5.44 (0.204) 


5.46 (0.006) 


5.52 (0.086) 






MRCI 


MP2 


4.80 (0.168) 


5.01 (0.003) 


5.13 (0.446) 


5.13 (0.004) 






Exp." 




4.5-4.6 


4.7-4.9 


5.8-6.1 






G 


TDDFT 


MP2 


4.84 (0.024) 


5.07 (0.142) 


5.22 (0.002) 


5.24 (0.014) 


5.42 (0.304) 






standard 


4.58 (0.001) 


5.04 (0.167) 


5.08 (0.002) 


5.25 (0.000) 


5.41 (0.313) 




MRCI 


MP2 


3.68 (0.001) 


4.74 (0.286) 


5.18 (0.002) 


5.21 (0.478) 


5.57 (0.014) 




Exp.*" 




4.4-4.6 


4.9-5.1 


5.5 


6.1-6.3 




T 


TDDFT 


MP2 


4.74 (0.000) 


5.22 (0.161) 


5.66 (0.000) 










standard 


4.74 (0.000) 


5.21 (0.156) 


5.63 (0.000) 








MRCI 


MP2 


4.63 (0.000) 


5.35 (0.434) 


5.74 (0.004) 








Exp." 




4.6-4.7 


5.6-6.1 


6.4 






C 


TDDFT 


MP2 


4.78 (0.049) 


4.84 (0.000) 


5.18 (0.002) 










standard 


4.70 (0.040) 


4.77 (0.000) 


5.07 (0.002) 








MRCI 


MP2 


4.69 (0.151) 


4.73 (0.002) 


5.69 (0.007) 








Exp."* 




4.5-4.6 


5.0-5.4 


5.6-6.1 







Average experimental excitation energies from Refs. "l48l|4i|50]'[5Tl|52l|5l[54]'=[50l|M[55l[5^ 



The transition densities between the ground and ex- 
cited states of the individual DNA bases were generated 
using TDDFT at the geometries of the bases optimized 
at MP2 level of theory, as described above. Before the 
actual calculation of the coupling elements could be car- 
ried out the transition densities and dipole moments ob- 
tained from ab initio calculations in an arbitrary coor- 
dinate system were transformed to the geometries of the 
bases in the studied DNA structures. This was carried 
out by defining the transformation superposing the plane 
defined by C6, Nl, and N3 atoms of pyrimidine or C6, 
Nl, and N9 atoms of purine bases in the arbitrary sys- 
tem with the plane defined by the corresponding three 
atoms of the base in the DNA structure. Subsequently, 
the transformation was applied to the three-dimensional 
grid holding the transition density and the dipole mo- 
ments. The quality of the fit as measured by the root- 
mean squared deviation between the atom coordinates of 
the two overlapped structures was very good. 



III. RESULTS AND DISCUSSION 

A. Individual Nucleobases 

The optimization of the standard nucleobase geome- 
tries [3S] obtained from X3DNA [33] at the MP2 level has 
very small effect on their geometries. The only notice- 
able difference is the out-of-planarity of the NH2 groups 
of adenine, guanine, and cytosine in the optimized ge- 
ometries. This is in agreement with previous theoretical 



studies of Shukla et al. ([36] and references therein) and 
experiment |37j where the amino groups of Ade, Gua, and 
Cyt were also found to be non- planar.. The root mean 
square deviations (RMSD) between the heavy atoms (ex- 
cluding hydrogens) of the original and optimized struc- 
tures is 0.045 A for adenine (Ade), 0.040 for guanine 
(Gua), 0.017 A for thymine (Thy), and 0.024 A for cyto- 
sine (Cyt). 



1. Excited- State Calculations 

The MP2 optimized structures of the 9H-purines and 
IH-pyrimidines were subsequently used to calculate the 
vertical excitation energies using time dependent density 
functional theory (TDDFT) at the PBEO/augTZVP level 
in gas phase. The results of the excited state calculations 
on Ade, Gua, Cyt, and Thy were also compared with the 
available experimental data and multireference configu- 
ration interaction (MRCI) calculations. 

Adenine. The lowest TDDFT calculated vertical sin- 
glet excitation energies of adenine, 5.00, 5.29, and 5.38 eV 
(Table |l]) correspond to the nn* and two closely spaced 
TTTT* transitions, respectively. While this order is in agree- 
ment with other DFT calculations (|36j). ab initio cal- 
culations at CASPT2 level ([38J) has shown the lowest 
excited state. Si, to be a tttt* state in agreement with 
the experiment. The UV-spectra of 9-methyladenine in 
stretched polymer poly(vinyl alcohol) films collected by 
Holmen (|39j) show two in-plane polarized transitions 
located at 4.55 and 4.81 eV. Contrary to the TDDFT 
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results the experimental data show the low-energy tttt* 
transition to carry less oscillator strength. The higher 
level ab initio calculations performed at MRCI level also 
predict the lowest energy state to be the light absorbing 
TTTT* state, calculated at lower energy, 4.81 eV, compared 
with TDDFT results. The second tttt* state is calcu- 
lated at 5.07 eV. In accordance with the experiment and 
CASPT2 data the MRCI predicts the higher energy tttt* 
transition to be more intense than the lower energy one. 
The most noticeable structural change, between the MP2 
optimized and standard geometry of Ade is the pyrami- 
dalization of the amine N in the former. The calculated 
transition energies at the TDDFT level imposed by the 
pyramidalization show a blue shift in the range of 0.10- 
0.15 eV for the flat structure. However, the separation 
between the two tttt* states and the character of the first 
excited state are not noticeable changed at this level of 
theory. 

Guanine. For the structure with the planar geome- 
try of the NH2 group (standard geometry) the computed 
lowest excitation energy at 4.59 eV is classified as the ira* 
transition. For this transition the configuration with the 
highest percentage weight, 99%, corresponds to HOMO 

LUMO, with the LUMO orbital being a a* locahzed 
at NH2 group. The lowest energy tttt* transitions for 
the flat structure are calculated at 5.04 and 5.41 eV. The 
pyramidalization of the NH2 group in the MP2 optimized 
structure causes the lowest energy transition to acquires 
some TT* character. It is now calculated at higher en- 
ergy 4.89 eV and defined by the configurations HOMO 
^ LUMO (91%) and HOMO -> LUMO+1 (5%). At this 
geometry the LUMO orbital is a mixture of tt* and cr* 
and the LUMO+1 is pure tt*. Similar mixing of tt* and 
cr* character was reported by Leszczynski and coworkers, 
who assigned corresponding transition for the nonplanar 
structure to the weak tttt* transition. In our calculations 
the two lowest energy tttt* transitions are calculated at 
5.07 and 5.42 eV and the transition at 4.89 eV is classified 
as TTa*. For this assignment of the tttt* transitions the 
difference in their calculated energies for the planar and 
pyramidal geometry of NH2 group is very small, below 
0.05 eV. A completely different situation was observed for 
the adenine, for which pyramidalization of the NH2 group 
caused blue shift of the tttt* transitions. The MRCI cal- 
culations yield the two lowest singlet vertical excitation 
energies at 4.24 and 4.34 eV which have mr* character. 
The lowest tttt* transitions at this level of theory are cal- 
culated at 4.76 and 5.24 eV. The latter tttt* transition 
has also larger calculated oscillator strength similarly to 
was was observed for adenine at this level of theory. 

Thymine. The ab initio and TDDFT calculations pre- 
dict the lowest energy nn* transition in vacuo for Thy in 
accordance with reported experimental data. The exci- 
tation energies calculated at the optimized and standard 
geometries of Thy are virtually the same for both the 
niT* and tttt* transitions. In an aprotic solvent thymine 
has the mr* state as the lowest singlet excited state [ID]. 
Present calculations at the TDDFT level, show dark mr* 



in: 



CaulsTiac mup ng bei|i."(^en 
'■V:ilfK> i-r;ri::k fini- 



\2oa 
laoD 



rvjo E:sckj9d arc o^ralkl 1hvriin&5 




5 a 

ni::l.'iri:v: 



■ - ■-■■IT. 



? 3 ; J 6 7 i 1.) 
r;;l;iii::i:;,'u'jK.;. 



FIG. 4: Coulombic coupling between the lowest tttt* transi- 
tion moments of Watson Crick AT base pair (left) and two 
stacked and parallel thymines (right) as a function of the dis- 
tance between the two bases. The distance is measured be- 
tween Nl and N3 atoms of A and T, respectively, for the AT 
pair (left) and between the centers of the mases of the two 
thymines (right). 



singlet excited state calculated at 4.74 eV, approximately 
0.5 below bright tttt* state calculated at 5.22 eV. At the 
MRCI level the relative order of this two transitions is 
the same and the calculated energies 4.63 and 5.35 eV 
of the mr* and tttt* transitions, respectively, are in good 
agreement with the corresponding TDDFT values. 

Cytosine. The TDDFT computed vertical singlet exci- 
tation energies of cytosine, shown in (Tablc|l]), predict the 
TTTT* state to be the lowest energy transition calculated at 
4.78 and 4.70 eV for MP2 and standard geometries, re- 
spectively. The non-planarity of the NH2 has only small 
effect on the energy of the tttt* transition inducing a red 
shift with a magnitude below 0.1 eV. The MRCI calcu- 
lations predicts the same order of the two lowest transi- 
tions, with TTTT* below mr* and their energies within 0.1 
eV of the corresponding TDDFT values (Table |l]). 



B. Coulombic Coupling 

The values of the Coulombic couplings between the 
lowest energy nn* transitions of the adenine and thymine 
and two 7r-stacked thymines as a function of distance 
between the bases (Fig. |4]) were calculated using the 
TDC and IDA methods. For the calculations the transi- 
tion densities and dipoles were those obtained from the 
TDDFT calculations on the MP2 optimized geometry of 
the basese. The comparison of the coupling elements ob- 
tained with the two methods, IDA and TDC, (Fig. ^ 
shows a good agreement at a separation between the 
bases larger than 5 and 6 A for the AT pair and two 
stacked thyminess, respectively. At a shorter separations, 
in the range of 3-4 A, which is typical for DNA struc- 
tures, the agreement between IDA and TDC is very poor 
with the differences between calculated couplings larger 
than 100% in case of AT pair. The aforementioned good 
agreement between IDA and TDC at a larger and poor 
agreement at shorter separations between nucleobases in- 
dicate that the shape and spatial extent of transition den- 
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sity (Fig. |4]) become important and cannot be neglected 
at distances between the bases typical for double helices 
DNA. The agreement between the two methods becomes 
very good in the limit of very large separation, (> 8 A). 

In order to compare the performance of the IDA and 
TDC method we calculated the coulombic couplings be- 
tween the lowest tttt* transitions of adenine and thymine 
bases for the interstrand Watson-Crick, intrastrand tt- 
stacked, and diagonal arrangement (Fig. [ij of the base 
pairs for the (^)io • {T)io oligomer in the idealized B- 
DNA geometry (Table ??) generated using the 3DNA 
program |33j using the IDA and TDC approximations. 
The corresponding band-structure as given by Eq. |6] is 
shown in Fig. [s] 

For the stacking and pairing distances correspond- 
ing to the idealized B-DNA geometry the coupling ele- 
ments calculated with the IDA approximation result with 
several-fold larger absolute values compared with the 
corresponding values calculated using the TDC method. 
The largest differences between the two methods are ob- 
tained for the couplings between the 7r-stacked adenines. 
For the idealized B-DNA geometry the coupling between 
two adenines located on the same strand calculated using 
IDA, 872 cm~^, is more than 5-fold larger compared with 
the value obtained using TDC 161 cm~^. The differences 
in the calculated couplings using the same two meth- 
ods for two stacked thymines are much smaller. For this 
base pair the Coulombic coupling calculated using IDA 
is equal to approximately 230 cm~^, more than twice the 
value of 101 cm~^ obtained with TDC. The absolute val- 
ues of the coupling elements between the second nearest 
neighbors located on the same strand are much smaller. 
At the IDA level of approximation the coupling between 



E(q)(cm \) 



44000 



42000 



40000 



3&000 



1 III 


- 








- 




ii=:^=^ 




Thy. band 

- 

1 


— ^ 



JT 
1 



JT 



the two adenines is only 57 cm"-'^ compared with 9 cm"-*^ 
calculated for the same base pairs using transition den- 
sity cubes. The coupling between the two thymine bases 
on the same strand is even smaller - approximately 3 and 
1 cm^^ for IDA and TDC methods, respectively. 

The values of the calculated coupling elements give 
hint on the relative exciton delocalization along the 
thymine and adenine strands of the (A)iq-(T)iq oligomer. 
The band structure in Fig. [5] shows that the mobility of 
the exciton along the thymine strand is low with both 
methods, IDA and TDC, giving resonable close results. 
The exciton mobility along the adenine strand, to the 
contrary, is quite different for the two methods. IDA, in 
this case, predicts more delocalized exciton states com- 
pared to TDA, which reflects larger magnitudes of the 
couplings calculated with the former method. 

The magnitudes of the couplings between bases lo- 
cated on different strands (Table ??), which belong to the 
Watson-Crick base pairs, are also larger when calculated 
with IDA method. The average Coulombic coupling for 
the Watson-Crick AT base pair calculated using IDA is 
230 cm^^ compared with 101 cm~^ obtained with TDC 
method. The magnitudes of the couplings between bases 
located on different strands, which does not belong to 
the Watson-Crick, the diagonal terms, are generally still 
quite large. Especially, the magnitude of the couphng 
between diagonal bases in the 5'-XY-5' direction (r"*") is 
comparable with the coupling for the Watson-Crick base- 
pairs (Table III). The values computed for the AT pair 
with the IDA and TDC methods are 163 and 109 cm"!, 
respectively. These values are almost twice and four 
times larger compared with the coupling between the cor- 
responding bases in the 3' direction (r~)calculated with 
the IDA and TDC methods, respectively. From the re- 
sults of Table ?? it is clear that the coupling elements are 
very sensitive with respect to the base sequence. The cal- 
culated Coulombic couphngs for the 7r-stacked arrange- 
ment of the adenine bases are by far the largest. Using 
IDA method the calculated coupling elements are more 
than two-fold larger than the corresponding values cal- 
culated with IDA for almost all arrangements. 

To investigate the effect of the structural fluctuations 
on the calculated couplings between the adenine and 
thymine basepairs we analyzed 4000 conformations of se- 
lected basepairs from the 80 ps molecular dynamics tra- 
jectories of the (A) 12 ■ (T)i2. The extracted 4000 snap- 
shots span the whole 80 ps simulations with each snap- 
shot taken every 20 fs. In Tabic III ) the average values of 



the couplings calculated using IDA and TDC methods for 
selected pairs of the {A) 12 ■ {T)i2 oligomer are listed. The 
comparison of these values with corresponding couplings 
calculated for basepairs in their idealized B-DNA geome- 
try (Table ??) shows some very interesting points. Com- 
paring the maximum and minimum values of the coupling 



FIG. 5: Band structure for AT B-DNA. Dashed curves: IDA, 
Solid: TDC. 



elements in Table III it can be seen that the magnitude of 
the couplings is very sensitive to the structural fluctua- 
tions observed in the MD simulations. The absolute value 
of the coupling can differ by as much as 1000 cm~^ for 
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TABLE II: Average values of the Coulombic couplings (cm"^) between the Adenine and Thymine bases of the (A)io-(T)io 
oligomer in the idealized BDNA geometry. The reported values are the averages for all corresponding base pairs with their 
standard deviations in parentheses. The hn, in, and correspond to coupling terms in Eqlsland Eq|6] 
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TT-stacked nucleobases. Other arrangements of the bases 
exhibit much smaller fluctuations of the couplings in the 
range of 300 cm~^ for Watson-Crick basepairs and still 
less for cross terms between bases on opposite strands 
and second nearest neighbors located on the same strand 
(Table III I . The fluctuations of the couplings observed for 
the intrastrand nearest neighbors show also more com- 
plex pattern compared to these for Watson-Crick base- 
pairs. An example is given by the couplings calculated 
with TDC method for the A9-A10 step (data not shown) 
of (A) 12 • (r)i2 oligomer. For this base pair the relatively 
larger and slower fluctuations on appoximately 20 ps time 
scale are superimposed on the rapid fluctuations. Sim- 
ilar slower fluctuations are observed for other 7r-stacked 
adenines but not Watson-Crick pairs. For both 7r-stacked 
adenines and thymines the Coulombic coupling seems to 
also depend on the location of a given pair along the 
chain. As can be inferred from data in (Table III) the 



base pairs closer to the 3' end, A9-A10 and T21-T22 show 
larger average values of the couplings compared with the 
corresponding values calculated for the base pairs closer 
to the 5' end. 

In spite of the large fluctuations of the couplings the 
values averaged over 4000 conformations are in good 
agreement with the corresponding couplings calculated 
for basepairs in standard B-DNA geometry. The best 
agreement is observed for the Watson-Crick base pairs 
while somewhat worse is seen for the 7r-stacked thymines. 
Interestingly, the average magnitudes of couplings for di- 
agonal base pairs in the 5' and 3' directions are very 
similar when calculated using IDA method. The average 
couplings of the A8-T16 and A9-T17 base pairs calcu- 
lated with this method for the 4000 conformations from 
the MD simulations are equal to 129 and 122 cm~^, re- 



spectively (Table III). The transition density cubes cal- 
culated couplings, however still show sensitivity the the 
direction with the average coupling in the 5' direction 
for A8-T16 , 77 cm"-'^, compared to only 33 cm~^ in 3' 
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FIG. 6: Eigenstate energy of (A)i2(T)i2 averaged over 240 
conformations from MD simulations obtained using TDC 
(filled circles) and IDA (filled squares). 



direction for A9-T17. 



C. Exciton states 

In this section we compare the properties of the excited 
states of the double-stranded DNA model (A) 12 -(T) 12 
calculated in the framework of the exciton theory where 
the off-diagonal elements of the exciton matrix were cal- 
culated using the ideal dipole approximation (IDA) and 
transition density cubes (TDC). The transition energies, 
oscillator strength, and the localization of the excited 
states were determined by diagonalization of the exci- 
ton matrix with the transition energies of the individual 
bases (diagonal terms) obtained from TDDFT calcula- 
tions using standard geometries of the nucleobases while 
the coupling elements (off-diagonal terms) were deter- 
mined using either IDA or TDC method. All properties 
were averaged over an ensamble of 240 conformations ex- 
tracted form the MD simulations. 

Fig. |6] shows the average values of the energies of the 
24 eigenstates of the (A) 12 • {T)i2 obtained with the cou- 
plings calculated using TDC and IDA method. Using the 
former method for calculations of the coupling elements 
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TABLE III; Average Coulombic coupling Vij (in cm~^) between selected bases of (A)i2-(T)i2. The values are averaged over 
4000 snapshots from an 80 ps MD simulation of DNA in water, cr is the r.m.s. deviation about the mean. Max. and Min. refer 
to the maximum and minimum value of the coupling over the entire simulation. 
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the energies of the lowest 12 eigenstates are clearly sep- 
arated from those of the remaining twelve eigenstates. 
The energy change between the two border eigenstates, 
<12> and <13>, amounts to almost 2500 cm"-'^ while 
the difference between the highest and the lowest energy 
eigenstates in each of the two sets is less than 1000 cm~^. 
The variations in the energy of a given eigenstate do not 
exceed 150 cm~^. The abrupt energy change between 



the two sets of eigenstates diminishes significantly when 
the IDA method is employed for calculations of the dipo- 
lar couplings (Fig.|6|. The difference between the average 
energies of the border eigenstates <12> and <13> is less 
than 700 cm"^ with variations in the energies of these two 
eigenstates equal to 420 and 350 cm~^, respectively. 

The spatial extent of the eigenstates was evaluated 
based on the inverse participation ratio {I PR) of a given 
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FIG. 7: Participation ratio of the 24 eigenstates of the 
(A)i2(T)i2. The values averaged over 240 conformations from 
MD simulations were calculated using TDC (solid lines) and 
IDA(dashed lines) methods. 
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FIG. 8: Contribution of the S2 states of thymine and adenine 
to the eigenstates of (A) 12 (T) 12 oligomer. 



eigenstate (Eq. |4| , which indicates the number of coher- 
ently bound chromophores [41]. The plot in Fig. |7] shows 
the average I PR values for each of the 24 eigenstates cal- 
culated using TDC and IDA to obtain the couplingele- 
ments of the Hamiltonian matrix. The average values of 
the inverse participation ratios for (A) 12 -(T) 12 obtained 
with TDC are in the range between 4.5 and 7.1. The 
corresponding values obtained using TDC method are 
slightly smaller for the first 12 eigenstates and slightly 
larger for the remaining 12 eigenstates and are in the 
range between 3.7 and 8.2. In both cases the IPR val- 
ues are much larger than one indicating delocalization 
of the excitation over several bases. Markovitsi et al. 
[H [131 E] showed that for the columnar aggregates of 
n identical chromophores, the maximum values of the 
normalized inverse participation ratio l/nLk is equal to 
0.7. Therefore, in case of the (A) 12 • {T)i2 oligomer a 
completely delocalized eigenstate over one strand of the 
double helix would have a participation ratio equal to 8.4. 
Contrary to what Bouvier reported for the (A) 20 -(T) 20 
and (dAdT)io-(dAdT)iQ oligomers, we found the partic- 
ipation ratios for all the eigenstates to be lower than 8.4. 
Therefore, a delocalization of the eigenstates over only 
adenosine or thymine chromophores but not both is ex- 
pected. As can be seen from the contribution of the S2 
transitions of adenine and thymine to the eigenstates of 
the (A)i2-(T)i2 (Fig. [s]) the lower energy eigenstates are 
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FIG. 9: Plot of the participation ratio of the eigenstate num- 
bers <13> and <22> as a function of energy determined for 
240 conformations of (A)i2(T)i2. 



localized almost completely on the S2 transition on the 
thymine, while the higher energy eigenstates are localized 
on the ^2 transition of adenosine. 

The inverse participation ratios of the eigenstates 
<13> and <22> as a funcion of energy are ploted in 
Figure [9j The IPR values for these two eigenstates of 
(A) 12 (T) 12 calculated for 240 conformations taken from 
the MD simulations show large fluctuactions in the range 
of 2 — 10. Despite the wide range of calculated IPR val- 
ues, as can be seen from the plots in Figure [9] the higher 
energy eigenstate, <22>, on average shows larger de- 
localization compared with the lower energy eigenstate, 
<13>, wheather the IDA or TDC method was used to 
calculate coupling elements. However, only for a hand- 
ful conformations the value of IPR exceeds the theoreti- 
cal value of 8.4 (indicated by a dashed line in Figure [9]) 
corresponding to the completely delocalized exciton over 
one strand of the (A) 12 (T) 12. This indicates that both 
eigenstates <13> and <22> which are localized on the 
transition associated with adenine remains localized on 
only one strand of the (A) 12 (T) 12 composed of adenine 
nucleobases. 

The average values of the oscillator strengths versus 
the average energies of the eigenstates for 240 confor- 
mations from MD simulation are plotted in Figure [TO] 
The total oscillator strength is distributed over a small 
number of eigenstates clustered in two bands. The 
first one comprise eigenstates <9> to <12> localized 
on thymine strand, and the second eigenstates <21> to 
<24> on adenine strand. The corresponding energies are 
around 39700 and 43000 cm for the off-diagonal cou- 
plings obtained using TDC and 40100 and 44000 cm"^ for 
off-diagonal terms calculated with IDA. These "bright" 
states correspond to the higher energy eigenstates built 
on the thymine and adenosine monomer transitions (Fig- 
ure 
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FIG. 10: Oscillator strength calculated for 24 eigenstates of 
the (A) 12 (T) 12. The off-diagonal elements of the exciton ma- 
trix were calculated using the TDC (top) and IDA (bottom) 
methods. The values were averaged over 240 conformations 
extracted from the MD simulations. 



ried by the eigenstates <10> and <22> and in both cases 
(IDA and TDA) the energies of these eigenstates are blue 
shifted with respect to transition energies of the individ- 
ual bases. The magnitude of the shift approximately 1000 
cm and less than 500 cm for IDA and TDC, respec- 
tively, reflects the differences in the couplings obtained 
with these two methods. 



IV. CONLUSIONS AND SUMMARY 

We have investigated the properties of the excited 
states of (A) 12 (T) 12 double helix calculated in the frame- 
work of the exciton theory. In our approach we combined 
the quantum mechanical calculations with the molecular 
dynamics simulations. The TDDFT calculations were 
employed to calculate the energies of the singlet excited 
states of the individual nucleobases. The transition mo- 
ments and densities of the Sq — > 15*2 transitions of adenine 
and thymine which correspond to the lowest energy tttt* 
transitions for these two bases were used to calculate the 
off-diagonal elements of the exciton matrix. The effect 
of the conformational changes were incorporated by av- 
eraging the calculated spectral properties of the double 
-stranded (^) 12 (7") 12 over large number of the conforma- 
tions extracted from the molecular dynamics simulations. 

The Coulombic couplings calculated using the IDA and 
TDC methods show a large deviations for the distances 
between chromophores typical for the DNA double he- 



lices. The magnitude of the couplings calculated with 
IDA being always larger than the corresponding values 
obtained with TDC. The agreement between the two 
methods is satisfactory only for the separations between 
the chromophores larger than 5 A. The largest difference 
between these two methods is observed for the 7r-stacked 
adenines in the standard B-DNA geometry for which the 
coupling calculated with IDA is over five times larger 
than the corresponding values calculated using TDC. The 
effect of the structural fluctuations on the calculated cou- 
pling elements is also very significant for both methods 
the values of the calculated coupling can change by an 
order of magnitude for different conformations of a given 
basepair. The difference between the smallest and largest 
coupling between the stacked adenines calculated using 
IDA for a given base pair can be as large as 1000 cm~^, 
smaller but still significant difference in the range of 300 
cm~^ was calculated using TDC. 

The properties of the excited states of the (A) 12(1") 12 
calculated in the framework of the exciton theory are af- 
fected to a different extent when the off-diagonal elements 
of the exciton matrix calculated using IDA and TDC 
methods. The eigenstates which carry the largest the 
oscillator strength, <10> and <22>, are slightly blue- 
shifted with respect to the transition energies of single 
nucleobases (Figure 10 1. The larger shift, approximately 
1000 cm"^, is observed for the exciton states obtained 
with the off-diagonal elements of the exciton matrix cal- 
culated using the IDA approximation, compared to less 
than 500 cm"^ obtained with the dipolar coupling cal- 
culated using TDC method. However, the delocalization 
properties of these eigenstates is similar for both IDA and 
TDC couplings. The IPR values of the "bright" eigcn- 
state <10> calculated with both IDA and TDC couplings 
are 5.5 and 6.0, respectively, while the corresponding IPR 
values of eigenstate <22> equal to 7.1 and 6.1 indicate 
only slightly large delocalization. Accordingly, compar- 
ing the IPRs obtained with TDC couplings the initial 
population of the bright eigenstates <10> and <22> by 
UV absorption will create exciton states which are de- 
localized over roughly 6 thymine and adenine bases, re- 
spectively. Upon relaxation the exciton states become 
more localized (Figures [9] and |7| as indicated by lower 
IPR values of the border eigenstates <1> and <13>. 
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